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Abstract 

We  analyzed  the  heartbeat  time  series  of  12  human  subjects  exposed  to  progressive  central 
hypovolemia  with  lower  body  negative  pressure.  Two  data  processing  techniques  based  on  wavelet 
transforms  were  used  to  determine  the  change  in  the  non-stationary  nature  of  the  time  series  with 
changing  negative  pressure.  Our  results  suggest  that  autonomic  neural  mechanisms  driving  cardiac 
interbeat  intervals  during  central  hypovolemia  go  through  various  levels  of  multifractility. 
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I.  INTRODUCTION 


In  the  past  two  decades  we  have  witnessed  an  explosion  in  the  biophysics  and  physiological 
literatnre  with  regard  to  the  identihcation  of  phenomena  having  long-term  memory  and 
probability  densities  that  extend  far  beyond  the  typical  tail  region  of  Gaussian  distributions. 
One  way  these  processes  have  been  classified  is  as  1/f-phenomena,  since  their  time  series 
have  spectra  that  are  inverse  power  law  in  frequency  or  their  probabilities  have  an  inverse 
power-law  distribution.  In  either  case  the  underlying  structure  is  fractal,  either  in  space, 


Q. 


time  or  both  ll||.  Physiological  signals,  such  as  the  heartbeat  interval  time  series  studied 
herein,  are  typically  generated  by  complex  self-regulatory  systems  that  process  inputs  with 
a  broad  range  of  characteristics.  Such  physiological  time  series  fluctuate  in  an  irregular 
and  complex  manner  such  as  shown  in  Fig.  1  and  the  statistics  of  these  fluctuations  often 


exhibit  self-affine  or  fractal  properties  l|.  Such  processes  do  not  posses_^  dominant  time 


scales,  but  they  can  frequently  be  characterized  by  fractal  dimensions 


Ivanov  et  al 


established  that  healthy  human  heartbeat  intervals,  rather  than  being  fractal,  exhibit 
multifractal  properties  and  they  uncovered  the  loss  of  multifractality  for  a  life-threatening 


condition  known  as  congestive  heart  failure.  Latka  et  al. 
cerebral  blood  flow  in  healthy  humans  is  also  multifractal. 


similarly  determined  that 


Recently,  Scafetta  et  al.  p  studied  how  the  fractal  and  multifractal  nature  of  the  human 
stride  interval  is  modified  by  changing  the  pace  velocity  in  both  free  and  metronomically 
constrained  walking  conditions  and  introduced  a  neuronal  super  central  pattern  generator 


(SCPG)  to  model  the  human  stride  interval  time  series  j^.  In  particular,  it  was  determined 
that  an  increase  of  neural  stress,  induced  by  a  faster  or  slower  pace  than  normal,  causes 
a  slight  increase  of  both  fractal  and  multifractal  properties  of  the  free  pace  stride  interval 
time  series  Isj.  In  the  present  analysis  we  used  lower  body  negative  pressure  (LBNP)  to 
induce  progressive  central  hypovolemia  and  determined  how  such  stimulation  of  autonomic 
neural  control  mechanisms  modified  the  statistical  nature  of  heart-beat  time  series.  A  LBNP 
is  applied  by  sealing  supine  subjects  at  the  level  of  the  iliac  crest  in  an  air-tight  chamber 
connected  to  a  vacuum  pressure  control  system.  The  magnitude  of  LBNP  is  carefully  regu¬ 
lated  to  induce  progressive  central  hypovolemia  and  autonomic  cardiovascular  compensatory 
mechanisms  are  evaluated.  LBNP  decreases  central  venous  pressure,  stroke  volume  and  car¬ 
diac  output  by  inducing  a  footward  fluid  shift,  and  increases  heart  rate  and  total  peripheral 
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resistance  through  parasympathetic  neural  withdrawal  and  sympathetic  neural  activation 
[^.  The  magnitude  of  changes  in  autonomic  neural  activity  can  be  estimated  non-invasively 
in  the  frmuency  domain,  but  common  analyses  incorporating  either  fast-Fourier  based  power 
spectra  js|  or  autoregressive  modeling  do  not  account  for  inherent  multifractality  or  non- 
stationarity  underlying  complex  cardiovascular  regulatory  mechanisms.  For  this  reason,  in 
the  present  study  we  applied  two  wavelet-based  techniques  to  study  the  fractal  and  mul¬ 
tifractal  dependence  of  the  cardiovascular  responses  to  progressive  central  hypovolemia  to 
determine  how  cardiovascular  stress  affects  fluctuations  of  the  data  on  a  hierarchy  of  time 
scales.  For  example,  the  analysis  allows  us  to  isolate  the  influence  of  the  breathing  cycle 


component  on 
gait  dynamics 


leart  rate.  We  determine  that  in  this  situation,  similar  to  the  situation  in 


30, 


an  induced  stress,  in  this  case  through  central  hypovolemia,  causes 
a  slight  increase  in  both  the  fractal  and  multifractal  properties  of  the  heartbeat  interval 
time  series.  We  suggest  that  such  changes  in  fractal  and  multifractal  properties  characterize 
neural  network  responses  of  biological  systems  to  physiological  stress. 

In  Sec.  II  we  present  the  heart  interbeat  interval  data  analyzed  herein;  Sec.  Ill  is  devoted 
to  a  review  of  the  techniques  of  analysis;  Sec.  IV  is  devoted  to  a  detailed  analysis  of  the 
heartbeat  interval  time  series  of  12  normal  subjects  and  in  Sec.  V  we  draw  some  conclusions. 


II.  HEART  INTERBEAT  DATA 

In  the  main  experiment  we  studied  twelve  healthy  men  between  the  ages  of  26  and  40 
after  they  had  been  informed  of  the  nature  of  the  experiment  and  had  signed  a  consent  form, 
approved  by  the  institutional  human  use  committee  of  the  U.S.  Army  Institute  of  Surgical 
Research,  Ft  Sam  Houston,  TX.  Subjects  were  positioned  supine  in  an  LBNP  chamber, 
sealed  at  the  waist  to  the  level  of  the  iliac  crest,  and  instrumented  with  a  standard  4- 
lead  electrocardiogram.  Because  respiratory  frequency  and  depth  may  seriously  confound 
interpretation  of  autonomic  rhythms  [l0|,  subjects  breathed  in  synchrony  with  a  metronome 
set  at  a  pace  of  15  breaths  per  minute  (0.25  Hz).  Electrocardiogram  data  were  sampled  at 
500  Hz  during  3-minute  periods  of  controlled  frequency  breathing  at  negative  pressure  of  0, 
-15,  -30,  -45,  -60,  and  again  at  0  mmHg.  Data  were  recorded  directly  to  computer  (Windaq, 
Dataq  Instruments,  Akron,  OH)  and  then  imported  into  a  customized  software  program  for 
analysis  [Cardio-Pulmonary  Research  Software  for  Windows,  (WinCPRS),  Absolute  Aliens 
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Ay,  Turku,  Finland].  R- waves  were  automatically  detected  from  continuous  waveforms  and 
marked  at  their  occurrences  in  time  to  obtain  heartbeat  periods  (tachogram)  as  shown  in 
Fig.  1. 

In  a  second  experiment,  conducted  at  the  Institute  of  Clinical  Physiology,  CNR-Pisa, 
we  compared  the  fluctuations  of  a  tachogram  with  the  contemporary  respiratory  signal  re¬ 
garding  one  person.  This  experiment  was  originally  done  with  thirteen  healthy  volunteers. 
These  volunteers  underwent  an  experimental  session  in  which  an  electrocardiographic  (ECG) 
derivation  and  the  respiratory  signal  were  simultaneously  and  continuously  recorded  at  1000 
Hz  sampling  frequency.  The  ECG  was  recorded  using  standard  electrodes,  while  the  respi¬ 
ratory  signal  was  detected  through  a  polymeric  piezoelectric  dc-coupled  transducer  inserted 
into  a  belt  wrapped  around  the  chest.  During  the  experimental  session  the  subjects  were 
comfortably  sitting  and  breathing  freely. 


III.  FRACTAL  AND  MULTIRESOLUTION  ANALYSIS  BASED  ON  WAVELETS 


The  analysis  we  apply  to  heartbeat  time  series  emphasizes  the  scaling  nature  of  the 
phenomenon.  In  the  next  two  subsections  we  briefly  review  the  two  techniques  we  use  in  our 
analysis.  The  hrst  method  determines  the  distribution  of  Holder  exponents  of  a  time  series, 
by  means  of  the  continuous  wavelet  transform  (CWT)  used  to  study 

the  fractal  and  multifractal  nature  of  a  time  series  and  has  already  been  used  to  analyze  the 
human  stride  interval  time  series  P|.  The  second  method  uses  the  wavelet  multiresolution 
analysis  (WMA)  of  a  time  series,  by  implementing  the  maximum  overlap  discrete  wavelet 
transform  (MODWT)  [1^.  The  latter  technique  provides  a  detailed  study  of  the  fluctuations 
of  a  time  series  on  a  scale-by-scale  basis. 


A.  Approximate  estimation  of  local  Holder  exponents  and  their  probability  dis¬ 
tribution. 


Mandelbrot  m  showed  that  many  natural  phenomena  are  described  by  self-affine,  cor¬ 
related,  scaling  time  series  and  as  we  mentioned  above,  heartbeat  data  fall  into  this  cate¬ 


gory.  The  sca^ 
noises  (fCn) 


ing  properties  of  the  fractal  noise  studied  by  Mandelbort,  fractional  Gaussian 


13j,  are  characterized  by  an  exponent  that  he  called  H  in  honor  of  Hurst. 
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Consequently,  if  X{t)  is  a  fractal  process  with  Hurst  exponent  H  and  c  is  a  constant,  then 
X  [t)  =  X{ct)/c^  is  another  fractal  process  with  the  same  statistics.  The  phenomenon  of 
fGn  has  a  spectrum  with  the  inverse  power-law  form 

P{f)  «  /-^  =  ^  ,  (1) 


where  /  is  the  frequency,  H  is  the  Hurst  exponent  and  Hq  is  the  average  of  the  Holder 
exponent  distribution  among  the  singularities  of  the  time  series.  We  can  also  express  a 
fractal  noise  through  the  equivalent  representation  of  the  autocorrelation  function  jl^ 


C[r)-  ocr 


(2) 


The  self-affine  property  expressed  by  m  and  the  relation  between  (3,  H  and  ho  are  theo¬ 
retically  valid  only  for  a  infinitely  long  monofractal  time  series.  Figure  2  shows  a  computer 
generated  realization  of  fGn  with  Hurst  exponent  FT  =  1  or  Holder  exponent  ho  ~  0,  also 
known  as  1//  noise  or  pink  noise.  This  type  of  noise  is  important  because  it  represents  a 
threshold  between  the  persistent-stationary  noise  (—0.5  <  ho  <  0)  and  the  non-stationary 
noise  (ho  >  0).  Moreover,  the  1//  noise  is  characterized  by  an  autocorrelation  function 
C(r)  that  does  not  depend  on  the  distance  r,  Eq.  (|2I).  The  random  noise  or  white  noise  is 
characterized  by  ho  ~  —0.5  whereas  the  random  walk  or  Brownian  motion  is  characterized 
by  ho  ~  0.5.  In  general,  the  higher  the  Holder  exponent,  the  smoother  the  time  series. 

As  seen  in  Figure  2,  a  fractal  noise  is  characterized  by  trends  and  discontinuities  that 
give  a  particular  geometric  shape  to  the  signal.  The  rapid  changes  in  the  time  series  are 


Q. 


called  singularities  of  the  signal  and  their  strength  is  measured  by  a  Holder  exponent 
Given  a  function  f{x)  with  a  singularity  at  xq,  the  Holder  exponent  h{xo)  at  such  a  point 
is  defined  as  the  supremum  of  all  exponents  h  that  fulfills  the  condition: 

\f{x)-Pn{x-xo)\<G\x-xo\^,  (3) 

where  Pn{x  —  xq)  is  a  polynomial  of  degree  n  <  h. 


The  Holder  exponent  of  a  singularity  can  be  evaluated  by  using  the  wavelet  transform 
|l^  3,Q|-  Wavelet  transforms  make  use  of  scaling  functions  that  have  the  property  of  being 
localized  in  both  time  and  frequency.  A  scaling  coefficient  s  characterizes  and  measures  the 
width  of  a  wavelet.  Given  a  signal  f{x),  the  continuous  wavelet  transform  (GWT)  of  f{x) 
is  defined  by 

OO 

"x  —  Xq 


if 


f{x)  dx  , 


(4) 
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where  the  kernel  'ip{u)  is  the  wavelet  hlter  centered  at  the  origin,  u  =  0,  with  nnit  width, 
s  =  1.  A  widely  adopted  choice  of  the  kernel  'ip{u)  is  the  so-called  Mexican  Hat,  that  is,  the 
second  derivative  of  a  Gaussian. 


The  wavelet  transform  can  be  used  to  determine  the  Holder  exponent  of  a  singularity  be¬ 
cause  the  wavelet  kernel  'ip{u)  can  be  chosen  in  such  a  way  as  to  be  orthogonal  to  polynomials 
up  to  degree  n,  that  is,  such  that  the  following  properties  are  fulhlled 


dx 


0  Vm,  0  <  m  <  n  . 


(5) 


In  fact,  if  (ini)  holds  true,  it  is  easy  to  prove  that  if  the  function  f{x)  fulhls  condition  (jH)),  its 
wavelet  transform  x  =  Xq  is  given  by 


-l-CXD 

m,r„(/)  =  f  il(u)  du  CC  Isl**™'  , 


(6) 


where  u  =  {x  —  Xq)/ s.  Therefore,  at  least  theoretically,  the  Holder  exponent  of  a  singularity 
that  is  localized  in  xq,  can  be  evaluated  as  the  scaling  exponent  of  the  wavelet  transform 
coefficient,  Ws^xoif)^  for  s  — >■  0. 

Even  if  Eq.  (jH}  can  be  evaluated  for  any  position  xn,  we  are  interested  only  in  the  cusp 
singularities  of  the  time  series.  Mallat  et  ah  lUl,  llS|  show  that  the  Holder  exponent 
of  these  singularities  can  be  evaluated  by  studying  the  scaling  exponent  h(xo)  along  the 
so-called  maxima  line  that  converges  towards  the  singularity.  The  maxima  lines  are  dehned 
by  the  extremes  of  the  wavelet  transform  coefficients  0  at  each  wavelet  scale  s.  Arneodo 
et.  al.  ji^.  1^  proved  that  WTMM  can  be  used  to  dehne  a  multifractal-like  formalism  that 
gives  the  stochastic  properties  of  the  singularities  of  a  fractal  or  multifractal  noise.  This 
methodology  has  been  recently  used  to  determine  the  multifractal  nature  of  many  signals. 


for  example,  that  for  human  heartbeats 


However,  the  above  method  presents  some  problems  of  stability  when  applied  to  observa¬ 
tional  data.  In  fact,  the  spectrum  can  be  corrupted  by  the  divergences  of  negative  moments 


21|  or  by  the  outliers,  that  is,  the  end  points  of  the  sample  singularities  11|.  Differ¬ 


ent  methods  have  been  suggested  to  remove  the  divergences  due  to  the  negative  moments 
of  the  multifractal  partition  function,  for  example,  by  chaining  the  wavelet  maxima  across 


scales 


Q  or. 


more  efficiently,  by  bounding  the  Holder  exponent  of  the  maxima  line  by  using 


the  slope  wavelet 


21|.  The  corruption  of  the  singularity  spectrum  due  to  the  outliers  is  more 


difficult  to  deal  with. 
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Struzik  [ijj  suggested  an  alternative  method  that  has  the  ability  to  determine  an  approx¬ 
imate  value  of  local  singularity  strength.  The  spectrum  may  then  be  evaluated  from  these 
approximate  values.  The  idea  is  to  estimate  the  mean  Holder  exponent  h  as  a  linear  £t  of 
the  following  equation 

log[M(s)]  =  h  log(s)  +  C  ,  (7) 

where  the  function  M(s)  is  obtained  via  the  partition  function,  and  the  mean  Holder  expo¬ 
nent  h  is  the  local  version  of  the  Hurst  exponent  H  □  ,  Mo.  p.c.el.  fa  a  .onfaaota, 


noise  with  Hurst  exponent  H,  we  have  h  =  H  —  1  because  the  Hurst  exponent  is  evaluated 
by  integrating  the  noise  [2,  ll^  •  Here  again  the  equality  is  only  rigorously  true  for  an  in- 
hnitely  long  monofractal  noise  data  set.  The  approximate  local  Holder  exponent  h{xo,  s)  at 
the  singularity  xq  can  now  be  evaluated  as  the  slope 

log(|hhs,xo(/)l)  -  (h  log(s)  +  C) 


h{xo,  s)  = 


(8) 


log(s)  -  log(sAr) 

where  sat  is  the  length  of  the  entire  wavelet  maxima  line  tree,  that  is,  the  maximum  available 
scale  that  coincides  with  the  sample  length  sn  =  N,  and  xq  belongs  to  the  set  f2(s)  of  all 
wavelet  maxima  at  the  scale  s  that  assume  the  value  Wg^xoif)- 


While  further  details  may  be  found  in  the  papers  by  Struzik 


ll|  and  in  Ref. 


Q. 


here  we 


focus  on  the  interpretation  of  the  final  output,  that  is,  the  distribution  of  estimated  local 


Holder  exponents  of  a  time  series  as  given  in  j^.  Fig.  3  shows  the  histogram  of  Holder 
exponents  obtained  using  a  computer-generated  data  set  with  a  Hurst  exponent  H  =  1 
shown  in  Fig.  2.  The  histogram  is  htted  with  a  normalized  Gaussian  function  of  the  type 


9{h)  = 


a 


exp 


{h  -  hpf 
2 


(9) 


that  captures  the  two  main  characteristics  of  the  time  series:  (a)  the  average  of  Holder 
exponent  ho  that  gives  an  estimate  of  the  fractal  nature  of  a  time  series;  (b)  the  width  of 
the  distribution  given  by  the  standard  deviation  a  that  measures  the  variability  of  the  local 
Holder  exponents,  that  is,  the  multifractality  of  a  time  series.  Finally,  as  Scafetta  et  al. 
^  explained,  a  monofractal  time  series  presents  a  particular  width  of  the  Holder  exponent 
distribution  that  depends  on  the  length  of  the  time  series,  as  Eq.  (jHl)  suggests.  Therefore, 
a  time  series  will  present  a  multifractal  nature  only  if  the  width  of  its  Holder  exponent 
distribution  is  larger  than  that  of  a  monofractal  time  series  of  the  same  length.  In  general, 
given  two  time  series  of  the  same  length,  if  ai  >  (J2,  we  may  say  that  the  time  series  ‘1’  has 
a  stronger  multifractal  character  than  time  series  ‘2’. 
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B.  Wavelet  multiresolution  analysis. 


The  second  wavelet-based  technique  that  we  implement  allows  a  detailed  study  of  the 


Qif 


fluctuations  of  a  time  series  on  a  scale-by-scale  temporal  basis.  The  MODWT  iia  is  almost 
independent  of  the  particular  family  of  wavelets  used  in  the  analysis  and  is  the  basic  tool 
needed  for  studying  the  multiresolution  analysis  of  time  series  of  N  data  points  via  wavelets. 


In  the  book  of  Percival  and  Walden, 


12j|.  the  reader  can  hnd  all  the  mathematical  details. 


In  this  paper  we  make  use  of  the  Daubechies  least  asymmetric  scaling  wavelet  filter  (LAS) 
that  looks  like  the  Mexican  Hat,  but  is  also  weakly  asymmetric;  a  fact  that  makes  LAS  filter 
more  malleable  than  the  Mexican  Hat.  The  WMA  via  MODWT  establishes  that  given  an 
integer  J  such  that  2"^  <  N,  where  N  is  the  number  of  data  points,  the  original  time  series 
represented  by  the  vector  X  can  be  decomposed  on  a  hierarchy  of  time  scales  represented 
by  a  smooth  part  plus  details  as  follows: 


^^  =  SJ  +  Y.D,, 

i=i 


(10) 


with  the  quantity  Sj  generated  by  the  recursion  relation 


Sj-i  —  Sj  -\-  Dj  .  (11) 

The  detail  Dj  of  Eq.  fllOj)  represents  changes  on  a  scale  of  r  =  2^,  while  the  smooth  Sj 
represents  the  smooth  wavelet  averages  on  a  scale  of  rj  =  2*^.  In  the  same  way  we  refer  to 
as  residuals  the  curves 

Rj  =  Y.r>,,  (12) 

i=i 

such  that  X  =  5*7  -I-  Rj. 

Fig.  4  shows  that  WMA  of  the  heartbeat  data  shown  in  Fig.  1  for  J=4.  The  detail 
curves  show  the  main  characteristic  fluctuation  that  characterize  each  time  scale  that  may 
be  expressed  in  the  number  of  beats,  as  done  in  Fig.  4,  or  in  physical  time  units.  We  shall 
use  these  detail  curves  to  establish  the  dependence  on  the  pressure  of  the  period  measured 
in  both  number  of  beats  and  in  physical  time  units  of  the  fluctuations  of  the  heartbeat  at 
different  wavelet  scale. 


IV.  ANALYSIS  OF  THE  HEARTBEAT  TIME  SERIES 


Fig.  1  shows  the  heartbeat  time  series  at  the  four  conditions  of  LBNP  for  a  typical 
individual.  In  this  work,  we  analyze  the  heartbeat  interval  time  series  for  12  individuals.  In 
the  time  series  data  shown  in  the  hgure  the  increase  in  the  average  heart  rate  with  increasing 
negative  pressure  is  clear,  as  is  the  return  to  normal  with  the  return  of  the  negative  pressure 
to  zero.  To  emphasize  this  effect,  in  Fig.  5  we  plot  the  average  heartbeat  interval  against  the 
pressure  for  12  individuals.  The  hgure  shows  that  the  average  heartbeat  interval  manifests 
a  steady  decrease  from  At  ~  980  ±  50  msec  at  P  =  0  mmHg  to  At  ~  710  ±  50  msec  at 
P  =  —60  mmHg.  At  the  recovery,  the  average  heartbeat  interval  is  slightly  larger  than  the 
initial  one;  we  measure  At  ~  1050  ±  50  msec.  The  error  bars  measure  the  average  standard 
deviation  of  the  variability  of  the  heartbeat  interval  for  each  individual  and  do  not  seem  to 
signihcantly  change  with  the  pressure. 

A  careful  visual  analysis  of  Fig.  1  suggests  that  by  increasing  the  negative  pressure, 
the  heartbeat  interval  time  series  becomes  smoother.  As  explained  in  the  previous  section, 
a  decrease  of  randomness  or  an  increase  of  the  smoothness  of  a  time  series  may  be  easily 
detected  by  the  Holder  exponents  of  the  time  series  that  measure  the  strength  of  the  singu¬ 
larities.  The  Holder  exponents  increase  with  the  smoothness  of  a  time  series.  Fig.  6  shows 
an  estimation  of  the  Holder  exponents  of  the  time  series  of  the  heartbeat  interval  time  series 
shown  in  Fig.  1.  We  use  the  algorithm  discussed  in  the  previous  section  to  obtain  these 
exponents.  Fig.  6  shows  that  the  Holder  exponents  signihcantly  increase  with  the  pressure, 
indicating  that  the  fractal  dimension  of  the  heartbeat  interval  monotonically  decreases  with 
negative  pressnre.  At  recovery,  the  fractal  properties  return  to  their  original  values. 

Fig.  7  depicts  the  average  Holder  exponent,  indicated  by  the  mean  value  ho,  and  by 
the  standard  deviation  a  for  the  12  individuals  against  the  pressure.  As  explained  in  the 
previous  section  the  standard  deviation  measures  the  variability  of  the  Holder  exponent 
values  that  can  be  related  to  the  multifractality  of  a  time  series.  The  data  indicates  that 
the  average  mean  value  has  a  steady  increase  from  ho  ~  0.177  ±  0.072  at  P  =  0  mmHg  to 
ho  ~  0.332  ±  0.104  at  P  =  —60  mmHg.  The  standard  deviation  is  indicated  by  the  error 
bars.  At  the  recovery  when  the  negative  pressure  returns  to  zero  we  get  ho  ~  0.157  ±  0.062. 
We  also  observe  that  by  increasing  the  negative  pressure  the  average  standard  deviation 
tends  to  increase.  The  values  we  measure  are  oq  =  0.072,  oys  =  0.062,  =  0.066, 


9 


=  0.077,  060  =  0.104  and  Oq  =  0.062  at  the  recovery.  As  anticipated,  an  increase  of  the 
standard  deviation  of  the  Holder  exponents  for  eqnal  length  time  series  may  be  interpreted 
as  an  increase  of  mnltifractal  properties  of  the  time  series.  Even  if  the  time  series  taken  for 
each  valne  of  the  pressnre  are  relatively  short  (approx.  200  data  points  for  P  =  0  mmHg 
and  approx.  300  data  points  for  P  =  —60  mmHg)  and,  therefore,  the  statistics  may  be 
not  optimnm.  However,  the  almost  monotonic  increase  of  the  standard  deviation  with  the 
pressnre  followed  by  an  abrnpt  decrease  at  the  recovery  may  reasonably  indicate  a  reliable 
increase  of  the  mnltifractality  of  the  heartbeat  interval  time  series  with  increasing  negative 
pressure. 

To  emphasize  the  properties  of  heartbeat  data  at  various  temporal  resolutions  we  need  to 
make  use  of  WMA  of  the  data.  In  Fig.  4  we  show  the  WMA  of  the  heartbeat  time  interval 
depicted  in  Fig.  1.  These  data  are  analyzed  using  the  MODWT  with  the  LAS  wavelet  hlter 
for  four  levels  of  resolution,  that  is,  until  the  wavelet  scale  J  =  4.  As  explained  in  the 
previous  section  the  smooth  curve  S4  captures  the  smooth  trend  of  the  fluctuation  at  the 
scale  of  r  =  2^  =  16  consecutive  heartbeats.  The  details  curves  D1-D4  capture  the  local 
fluctuations  of  the  signal  at  each  of  the  four  time  scale  r  =  2*.  This  means,  for  example, 
that  D1  captures  the  fluctuations  of  the  original  signal  with  a  period  that  is  approximately 
between  2  and  4  beats.  We  recall  that  the  shortest  possible  period  that  can  be  detected  in 
a  time  series  is  2  time  units.  The  detail  D2  captures  the  fluctuations  with  a  period  in  the 
interval  between  4  and  8  beats,  and  similarly  for  the  following  details.  According  to  Eq. 
(uni)  the  sum  of  the  four  details  D1-D4  and  the  smooth  curve  S4  gives  the  original  datasets 
shown  in  Fig.  1. 

Fig.  4  stresses  the  great  utility  of  analyzing  times  series  of  complex  systems  through  the 
WMA  via  the  MODWT.  In  fact,  the  wavelet  sensibility  to  the  local  changes  of  the  signal 
allows  us  to  easily  extract  information  that  is  impossible  to  obtain  by  using,  for  example, 
the  Fourier  transform  since  the  latter  averages  the  changes  of  the  entire  signal  at  each 
Fourier  transform  frequency.  For  example.  Fig.  4  clearly  shows  that  the  amplitude  of  the 
local  decomposition  is  not  constant  but  fluctuates  in  time.  In  particular  we  notice  that  the 
amplitude  of  the  fluctuations  captured  by  the  details  D3  and  D4  does  not  seem  to  change 
systematically  with  the  pressure.  Therefore,  we  can  conclude  that  the  change  of  pressure  do 
not  influence  signihcantly  the  slow  fluctuations  with  a  period  larger  than  8  heartbeats. 

On  the  contrary,  the  details  D1  and  D2  clearly  show  a  signihcant  reduction  of  the  ampli- 
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tude  of  the  fast  fluctuation  with  increasing  negative  pressure.  At  the  recovery,  the  amplitude 
of  the  fluctuations  returns  to  a  value  similar  to  the  initial  one  for  P  =  0  mmHg.  The  re¬ 
duction  of  the  amplitude  of  the  fast  fluctuations  manifests  itself  with  the  increase  of  the 
smoothness  of  the  heartbeat  interval  time  series  that  we  have  already  noticed  to  occur  by 
increasing  the  negative  pressure,  as  Fig  1  shows.  We  have  also  noticed  that  an  increase  of 
the  smoothness  yields  an  average  increase  of  the  local  Holder  exponents,  that  are  associated 
with  the  fractal  dimension  of  a  time  series,  as  explained  earlier  and  shown  in  Figs.  6  and  7. 
Fig.  4  shows  the  importance  of  using  WMA  because  we  are  able  to  isolate  the  component 
of  the  signal  of  the  heartbeat  interval  time  series  that  presents  dependence  on  the  pressure 
and  identihes  it  with  the  changes  captured  by  the  wavelet  details  D1  and  D2. 

Because  each  detail  curve  captures  the  variation  of  the  signal  inside  a  frequency  or  period 
range,  it  is  possible  to  also  study  how  the  average  period  of  the  fluctuations  at  a  particular 
wavelet  scale  are  affected  by  the  pressure.  Fig.  8  shows  the  average  mean  period  measured 
in  number  of  beats  of  the  fluctuation  associated  to  each  detail  curve  obtained  by  the  WMA 
for  12  subjects.  We  stress  that  these  calculations  are  made  without  any  kind  of  interpolation 
between  the  data.  So  the  natural  unit  of  the  time  sequence  is  the  number  of  heartbeats. 
To  get  the  results  in  physical  time  units  we  have  to  convert  the  results  obtained  with  the 
heartbeat  number  in  physical  time  units  by  using  the  original  R-R  sequence  or  we  may 
attempt  a  direct  calculation  of  WMA  after  an  interpolation  procedure  that  generates  an 
even  time  series  in  a  physical  units.  However,  interpolation  alters  data,  so  we  prefer  the  hrst 
procedure.  Fig.  9  shows  the  average  mean  period  measured  in  physical  time  units  of  the 
same  fluctuations.  We  decided  to  plot  the  results  concerning  Dl,  R2,  D3  and  D4  curves. 
We  stress  that  we  are  considering  the  residual  R2  that  is  given  by  R2  =  Di  +  D2  according 
to  m  rather  than  the  D2  curve.  We  prefer  to  use  the  residual  R2  instead  of  the  detail 
D2  because  the  range  from  2  to  8  heartbeats  is  supposed  to  capture  the  high  frequency 
component  of  the  R-R  signal  that  is  coupled  with  the  respiratory  signal.  The  residual  R2 
fully  contains  such  information  that,  instead,  is  split  between  Dl  and  D2  curves  with  Dl 
being  a  kind  of  harmonic  of  D2.  For  this  experiment  we  recall  that  the  subjects  breathed  in 
synchrony  with  a  metronome  sets  a  pace  of  15  breaths  per  minute  (0.25  Hz).  So  we  expect 
to  hnd  such  a  frequency  in  the  R2  curves.  The  detail  D3  captures  the  range  from  8  to  16 
heartbeats  and  corresponds  to  the  low  frequency  component  of  the  R-R  signal.  Higher  levels 
of  details,  that  are  included  in  the  S3  smooth  curve,  capture  the  very-low  component  of  the 
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interbeat  signal. 

The  results  plotted  in  Figs.  8  and  9  concerning  the  detail  curves  D3  and  D4  show  a 
more  erratic  behavior  with  large  errorbars  and  their  interpretation  remains  unclear.  Fig.  8 
shows  an  increase  of  the  heartbeat  number  for  cycles  captured  by  D3  curve  with  increasing 
negative  pressure.  Fig.  9  shows  that  in  physical  time  unit  the  period  of  these  D3  cycles 
remains  almost  constant  at  almost  10  seconds,  that  is,  at  a  frequency  of  0.1  Hz,  while 
the  period  of  the  cycles  associated  with  D4  detail  curve  tends  to  decrease  with  increasing 
negative  pressure. 

The  results  concerning  the  fast  fluctuations,  curves  D1  and  R2,  show  more  regular  pat¬ 
terns.  In  particular  we  notice  that  by  increasing  the  negative  pressure.  Fig.  8  shows  that  the 
average  mean  period  of  the  fluctuations  measured  in  heartbeat  number  is  almost  constant 
for  details  Dl,  D3  and  D4.  Instead,  the  results  concerning  the  R2  residual  shows  a  sensitive 
steady  increase  of  the  average  number  of  heartbeats  from  almost  4.1  heartbeats  per  cycle  to 
4.9  heartbeats  per  cycle  with  increasing  negative  pressure  from  P  =  0  mmHg  to  P  =  —60 
mmHg.  At  recovery  the  average  number  of  heartbeats  returns  to  4.1  heartbeat  per  cycle.  By 
looking  at  the  same  results  expressed  in  physical  time  units,  see  Fig.  9,  we  observe  that  the 
average  period  of  the  fluctuations  captured  by  the  R2  residual  curve  is  almost  constant  with 
the  pressure.  We  measure  4  seconds  that  correspond  to  a  frequency  of  0.25  Hz.  Instead,  the 
average  period  of  the  fluctuations  captured  by  the  Dl  detail  curve  shows  a  steady  decrease 
from  3.2  to  almost  2.1  seconds  with  increasing  negative  pressure.  The  different  behaviors 
characterizing  the  wavelet  detail  Dl  and  R2  depicted  in  Figs.  8  and  9,  suggest  that  the 
two  curves  can  be  partially  decoupled  because  they  respond  differently  to  the  increases  in 
negative  pressure. 

Fig.  10  shows  a  comparison  between  the  R2  residual  curve  of  the  heartbeat  interval  and 
respiratory  signals  of  a  person  for  60  seconds  under  normal  conditions.  These  data  were  ob¬ 
tained  in  a  second  experiment.  Figure  10  shows  a  very  strong  correspondence  between  the 
fluctuations  captured  by  the  R2  residual  curve  and  the  respiratory  cycle.  This  comparison 
suggests  that  the  R2  residual  curve  has  the  physiological  meaning  of  describing  the  fluctua¬ 
tions  of  respiratory  sinus  arrhythmia.  In  fact,  as  Fig.  9  shows,  the  fluctuations  captured  by 
the  R2  residual  have  the  frequency  of  0.25  Hz  of  the  respiratory  signal  that  is  kept  constant 
during  the  hypovolemia  experiment. 
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V.  SUMMARY  AND  CONCLUSIONS 


The  analysis  of  the  data  suggests  that  progressive  central  hypovolemia  induced  with 
LBNP: 

1)  Causes  decreases  of  R-R  intervals  and,  therefore,  increases  of  heartbeat  rate. 

2)  Decreases  heart  rate  variability  and  therefore  smoothes  out  the  heartbeat  time  series. 

3)  Reveals  that  two  measures  of  the  decrease  in  the  heartbeat  rate  variability  with  in¬ 
creasing  negative  pressure  are  an  increase  in  the  Holder  exponent  and  a  decrease  in  the 
fractal  dimension  of  the  time  series. 

4)  Reveals  that  accompanying  the  increase  in  the  Holder  exponent,  with  increasing  neg¬ 
ative  pressure,  is  an  increase  in  the  multifractality  of  the  time  series  as  measured  by  the 
width  of  the  Holder  exponent  distribution.  This  non-stationary  aspect  of  the  data  was  also 
observed  in  the  study  of  human  gait  data,  when  the  walkers  were  put  under  physical  stress 


by  asking  them  to  walk  faster  or  slower  than  normal 

5)  Reveals  that  independent  measures  of  the  heartbeat  rate  response  to  the  increases  in 
negative  pressure  are  the  wavelet  transform  details  D1  and  D2.  Both  these  non-periodic 
functions  display  a  monotonic  decrease  in  amplitude  with  increasing  negative  pressure. 

6)  Suggests  that  the  amplitude  of  the  higher-order  details,  D3  and  D4,  are  independent 
of  the  changes  in  the  negative  pressure.  However,  the  period  of  the  cycles  captured  by  D4 
tends  to  decrease  with  increasing  negative  pressure.  This  suggests  that  part  of  the  spectral 
energy  of  the  very-low  frequency  (/  <  0.04  Hz)  component  of  the  RR  signal  moves  toward 
the  low  frequency  (0.04  <  /  <  0.15  Hz)  component  of  the  same  signal. 

7)  Suggests  that  the  wavelet  transform  residual  curve  R2  is  driven  by  the  respiratory  cycle. 
The  residual  R2  captures  the  high  frequency  component  of  the  RR  signal  (  0.15  <  /  <  0.4 
Hz) .  This  is  an  indication  that  the  wavelet  transform  is  able  to  sort  through  the  fluctuations 
of  the  cardiac  time  series  and  highlight  particular  physiological  phenomena. 

These  separate  observations  made  from  the  analysis  of  the  heartbeat  time  series  stressed 
by  central  hypovolemia  are  interesting  in  themselves,  but  they  suggest  a  much  more 
signihcant  hnding.  They  suggest  that  the  influence  of  stress  on  a  physiological  system 
is  manifest  in  the  non-stationary  character  of  the  time  series  and  that  the  degree  of 
non-stationarity  can  be  measured  through  the  level  of  fractality  and  multifractality  of  the 
time  series  regarding  the  stressed  system. 
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FIG.  1:  Typical  response  of  the  heartbeat  period  at  pressure  change. 
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FIG.  2:  Gomputer-generated  fractal  noise  (FGN)  with  Hurst  coefficient  H  =  also  known  as  1// 
noise  or  pink  noise  is  shown. 
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FIG.  3:  Histogram  and  probability  density  estimation  of  the  Holder  exponents  for  the  computer¬ 
generated  fractal  noise  with  Hurst  coefficient  H  =  1  shown  in  Fig.  2.  The  fitting  curve  is  a 
Gaussian  0  centered  in  /iq  =  0.007  ±  0.001  and  with  a  width  a  =  0.051  ±  0.001.  The  mean 
Holder  exponent  given  by  Eq.  0  is  h  =  —0.004  ±  0.008  and  correspond  to  the  maximum  of  the 
distribution. 
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Wavelet  multiresolution  analysis  for  J  =  4  for  the  data  shown  in  Fig.  1. 
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FIG.  5:  Average  heartbeat  period  as  function  of  the  pressure  for  the  12  datasets. 
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FIG.  6:  Local  Holder  exponent  for  the  singularities  of  the  hearbeat  interval  time  series  shown  in 
Fig.  1. 
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FIG.  7:  Holder  exponent  distribution  for  the  heartbeat  period  time  series  as  function  of  the  pres¬ 
sure.  The  error  bars  indicate  the  average  standard  deviation  of  the  Holder  exponent  distribution 
for  the  12  datasets. 
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FIG.  8:  Average  mean  period  measured  in  number  of  heartbeats  of  the  fluctuation  associated  to 
Dl,  R2,  D3,  and  D4  detail  curve  obtained  by  the  WMA  for  J  =  4  for  the  12  datasets. 
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FIG.  9:  Average  mean  period  measured  in  physical  time  units  of  the  fluctuation  associated  to  Dl, 
R2,  D3,  and  D4  curve  obtained  by  the  WMA  for  J  =  4  for  the  12  datasets.  Note  that  the  period 
associated  to  the  residual  R2  remains  4  second  under  different  intensity  of  negative  pressure.  Four 
seconds  is  the  period  of  the  breathing  cycle. 
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FIG.  10:  Comparison  between  R2  detail  curve  of  the  heartbeat  interval  signal  and  the  contempo¬ 
rary  respiratory  signal  of  a  subject. 
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